;+
; NAME:
;       ASSOC_BGPS_BOLOCAT
;
; PURPOSE:
;       Creates the IDL save file local/BGPS_v102_map_locations.sav
;       which associates each Bolocat object with a position in the
;       BGPS release maps.  Selects BGPS mosaic based on label maps
;       produced by Bolocat.
;
; CATEGORY:
;       distance-omnibus utility
;
; CALLING SEQUENCE:
;       ASSOC_BGPS_BOLOCAT [,START=start]
;
; INPUTS:
;       NONE
;
; OPTIONAL INPUTS:
;       START -- Bolocat object number to start with [Default: first].
;
; KEYWORD PARAMETERS:
;       NONE
;
; OUTPUTS:
;       NONE (Creates the IDL save file local/BGPS_v102_map_locations.sav)
;
; OPTIONAL OUTPUTS:
;       NONE
;
; MODIFICATION HISTORY:
;
;       Created:  08/25/10, TPEB -- Initial (undocumented) version.
;       Modified: 07/14/11, TPEB -- Added documentation, code cleanup.
;       Modified: 07/19/11, TPEB -- Added check for correct PPBEAM
;                                   FITS keyword in the BGPS header.
;       Modified: 01/12/12, TPEB -- Additional comments and code
;                                   cleanup.
;       Modified: 07/23/12, TPEB -- Added smoothed BGPS map name to
;                                   the structure.
;
;-

PRO ASSOC_BGPS_BOLOCAT, START=start
  
  COMPILE_OPT OBSOLETE
  
  IF n_elements(start) EQ 0 THEN start = 0L ELSE start = long(start)
  
  ;; Read in Bolocat v1.0.1
  s = read_bgps_csv('bgps_distance_database_blank.csv',csv,/VERBOSE)
  
  ;; Find all the BGPS v1.0.2 data and label maps
  SMOOTH_DIR = './local/bgps/other/'
  LABEL_DIR  = './local/bgps/label/'
  MAP_DIR    = './local/bgps/maps/'
  bgps_maps  = FILE_SEARCH(MAP_DIR  +'*.fits', COUNT=n_map)
  bgps_list  = FILE_SEARCH(LABEL_DIR+'*.fits.gz', COUNT=n_list)
  bgps_smoo  = FILE_SEARCH(SMOOTH_DIR+'*_13pca_smoothmap50.fits', COUNT=n_smoo)
  
  message,string(n_map,n_list,n_smoo,format=$
                 "('N_MAP: ',I0,'  N_LIST: ',I0,'  N_SMOO: ',I0)"),/inf
  
  IF (n_list NE n_map) || (n_list NE n_smoo) || (n_map NE n_smoo) THEN $
     message,'Mismatch between labels, maps, and smoothed images!'
  
  ;; Make structure for holding image association locations
  bgps = replicate( CREATE_STRUCT('cnum',long(0),'filename','','labelname','',$
                                  'smooname','','xpos',long(0),'ypos',long(0),$
                                  'glon',double(0),'glat',double(0)), csv.nrows)
  
  ;; Fill in various structure elements
  bgps.cnum = s.cnum
  bgps.glon = s.glon_peak
  bgps.glat = s.glat_peak
  q = '"'
  
  ;; Read the label images & headers to memory (for faster processing)
  FOR j=0L, n_list-1 DO BEGIN
     
     command = string(j,q,bgps_list[j],q,format=$
                      "('image',I0,' = readfits(',A0,A0,A0,',hdr)')")
     errcode = Execute(command)
     
     command = string(j,format="('extast,hdr,astr',I0)")
     errcode = Execute(command)
     
     ;; Also, check the FITS header for the correct value of PPBEAM,
     ;;   and write the corrected value, if necessary.
     IF check_bgps_ppbeam( hdr ) THEN BEGIN
        command = string(q, bgps_list[j],q,j,format=$
                         "('writefits, ',A0,A0,A0,',image',I0,',hdr')")
        errcode = Execute(command)
     ENDIF
     
  ENDFOR
  
  ;; Loop through each BGPS source
  FOR i=start, csv.nrows-1 DO BEGIN
     
     count = 0
     
     ;; Loop through each label image
     FOR j=0L, n_list-1 DO BEGIN
        
        command = string(j,format=$
                         "('ind = WHERE(image',I0,' EQ s[i].cnum, nind )')")
        errcode = Execute(command)
        
        ;; When we find a match...
        IF (nind NE 0) THEN BEGIN
           count += 1
           
           IF ((i+1) MOD 100) EQ 0 THEN $
              print,s[i].glon_peak,s[i].cnum,$
                    nind,'   ',string(j,format="('Image',I0)"),'  ',bgps_maps[j]
           
           ;;=======================================================
           ;; Get the BGPS image filename, and place in structure
           ;; Get the peak position of the BGPS source in the image
           command = $
              string(j,format=$
                     "('ad2xy,s[i].glon_peak,s[i].glat_peak,astr',I0,',x,y')")
           errcode = Execute(command)
           
           command = string(j,format="('imsize = size(image',I0,', /DIM)')")
           errcode = Execute(command)
           
           x = long(round(x))
           y = long(round(y))
           
           IF( x LT 0  OR x GE imsize[0] ) THEN BEGIN
              print,string(bgps[i].cnum,x,format=$
                           "('BGPS #',I4,'  WARNING!!!  XPOS: ',I0)")
              count -= 1
           ENDIF ELSE BEGIN
              bgps[i].filename  = $
                 strmid(bgps_maps[j],strpos(bgps_maps[j],'/',/REVERSE_SEARCH)+1)
              bgps[i].labelname = $
                 strmid(bgps_list[j],strpos(bgps_list[j],'/',/REVERSE_SEARCH)+1)
              bgps[i].smooname = $
                 strmid(bgps_smoo[j],strpos(bgps_smoo[j],'/',/REVERSE_SEARCH)+1)
              bgps[i].xpos = x
              bgps[i].ypos = y
           ENDELSE
        ENDIF
        
     ENDFOR
     
     IF count NE 1 THEN $
        message,string(bgps[i].cnum,count,format=$
                       "('BGPS #',I4,' has ',I0,' matches!')"),/inf
     
  ENDFOR  
  
  save,bgps,filename='local/BGPS_v102_map_locations.sav',/verbose
  
  ;; Clean up the memory
  undefine,bgps
  
END
